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Emergent behaviors occur in a vast array of systems across many scales, and are of fundamental physical 
importance because of the intrinsic difficulty in linking microscopic system properties to macroscopic behav¬ 
iors. Here we study the emergent self-assembly behavior of model systems of recently synthesized families 
of concave dimpled hard spheres, or lock-and-key colloids. We find that as dimple size increases each family 
exhibits a crossover from a structure that does not reflect the particle symmetry to one that does and, surpris¬ 
ingly, the point at which this crossover occurs is approximately independent of the particle symmetry. Using 
a combination of numerical and analytic techniques we study systems at infinite and finite pressure, and find 
different common control parameters in each limit. Our results suggest there exists a set of experimentally re¬ 
alizable colloidal systems that exhibit complex emergent behaviors that can be traced to a common underlying 
microscopic control parameter. 


I. INTRODUCTION 


Physical systems from flocking animals [fTl to heavy 
fermion materials m exhibit emergent macroscopic behav¬ 
iors that are intrinsically difficult to predict from their micro¬ 
scopic properties El. Colloidal systems are fertile ground for 
investigating emergent behaviors for three reasons: a wide 
range of colloidal systems manifest emergent behaviors be¬ 
cause they have entropy-driven phase behavior 110; there 
is a broad range of techniques for manipulating colloids ex¬ 
perimentally (Ml; and colloidal systems can be investigated 
with a variety of techniques including in situ confocal mi¬ 
croscopy cni. Here, we study how microscopic properties 
affect the emergent macroscopic behavior of a model family 
of anisotropically shaped colloidal particles. 

Anisotropic colloids that exhibit emergent directional en- 
tropic forces (DEEs) o when crowded in a thermal sys¬ 
tem are entropically patchy particles. ifTll Entropic patches are 
shape features that promote local packing into target motifs 
and drive systems to order globally at high density. mi Previ¬ 
ously, we found that global ordering across a variety of sys¬ 
tems was consistent with entropic patch strengths on the or¬ 
der of a few k^T. EdD However, if we desire particles with 
entropic patches of this strength, we do not currently have a 
complete set of rules for understanding how the particle geom¬ 
etry needs to be controlled to get them. da Due to the emer¬ 
gent nature of the underlying mechanism, establishing quanti¬ 
tative control is difficult. Bimiia Because of this difficulty, 
the rules we previously proposed for controlling DEEs were 
heuristic. mi Here, we attempt to make the rules for control¬ 
ling emergent DEEs more quantitative for a speciflc case of 
families (see Eig.[T]) of anisotropic lock-and-key colloids lfT3]| . 

Lock-and-key colloids are a convenient model system for 
investigating how to use geometry to control DEEs. More 
importantly, a number of theoretical EHIll and experimen¬ 
tal investigations iniiiMa have shown that lock-and-key 
colloids can form interesting structures, and recently synthe¬ 
sized families of multiply dimpled particles |[T9l [24l have not 
been previously simulated. Here, we use computer simulation 




0.00 0.25 


0.50 0.75 


4 . 

1.00 



9 IB 


P2 P3 P4t 


P4s P6 


FIG. 1. (a) Change in shape of a tetrahedral tetravalent particle from 
convex to concave as dimpling amount is increased from 0 to 1. (b) 
Parameterization of dimpled concave spheroidal particles, (c) Biva¬ 
lent, trivalent, tetrahedral tetravalent, square tetravalent and hexava- 
lent model spheroidal particles. 


and analytic calculations to study both packing (inflnite pres¬ 
sure) and assembly (flnite pressure) behaviors of families of 
multiply-dimpled particles that differ by dimple number and 
symmetry as a function of dimple size. We And that in both 
cases there are common, approximate microscopic parame¬ 
ters, based on particle geometry, that correlate with changes 
in bulk system behavior. However, the microscopic parame¬ 
ters are different in the dense packing and assembly limits. In 
the dense packing limit, a combination of analytic and numer¬ 
ical calculations show that (i) as dimple size increases, there 
is a crossover from body-centered-cubic (BCC) or -tetragonal 
(BCT) phases to phases that reflect the number of dimples and 
particle symmetry, and surprisingly (ii) this point of crossover 
occurs approximately independent of particle symmetry and 
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number of dimples. Similarly, in the assembly limit, we find 
that as dimple size increases in each family there is (i) a 
crossover from the face-centered-cubic (FCC) structures that 
are characteristic of hard sphere assembly fTf\ to BCC struc¬ 
tures that are characteristic of soft sphere assembly (281 . fol¬ 
lowed by (ii) a second crossover from BCC to behavior that is 
determined by the particle symmetry and number of dimples. 
Surprisingly (iii) this second point of crossover also depends 
on a quantity that is approximately independent of particle 
symmetry, though it is different from the control parameter 
in the dense packing limit. We find that in the dense packing 
limit, the crossover in behavior is determined by a quantity 
that depends simultaneously on the geometry of all of the par¬ 
ticle features, whereas in the assembly limit it is determined 
by the volume of individual dimples. We interpret these re¬ 
sults in light of previous work 13 arguing that shape entropy 
maximization causes entropically driven colloidal systems to 
optimize local packing. Our results provide detailed quantita¬ 
tive evidence on how shape features control the macroscopic 
phase behavior of colloidal systems, while providing a con¬ 
crete example of a set of experimentally realizable systems 
exhibiting changes in emergent macroscopic behaviors that 
are correlated with changes in a common microscopic control 
parameter. 


II. MODEL AND METHODS 

A. Model Systems 

Particles, shown schematically in Fig. [Jc) and denoted by 
Fhd, are made of a central sphere P that is symmetrically dim¬ 
pled by subtracting no valence spheres, all of the same radius 
as the central sphere (hereafter taken to be one). The valence 
sphere positions are chosen to be related by discrete symme¬ 
tries, and are equidistant from the central sphere as shown in 
Fig. [Jb). The use of equal radii suggests a connection be¬ 
tween the dimple volume and free volume gained by the rest 
of the system when particles bind entropically (S) [TTl due to 
shape complementarity. We study five families of particles 
that satisfy two criteria: they are multi-dimpled (so that we 
can differentiate between effects that are driven by individ¬ 
ual shape features and effects that are driven by overall parti¬ 
cle shape) and have discrete symmetries that are commensu¬ 
rate with some lattice symmetry group. We denote the parti¬ 
cles as: bivalent (P2), trivalent (P3), tetrahedrally tetravalent 
(P4t), square tetravalent (P4s), and hexavalent (P6). For each 
of these concave particles, there is a one-parameter family of 
shapes determined by the distance d between the central and 
valence spheres which is maximum when the valence spheres 
are tangent to the central sphere and minimum when the va¬ 
lence spheres are tangent to each other. 

Dimples are pairs of spherical caps bounded by the intersec¬ 
tion of the central and valence spheres, and have the following 


TABLE I. Particles from the Concave Spheres Family 


Particle Description Arrangement 

Domain 

Dimpling Amount 

P2 

Bivalent 

Linear 


(A-d^) 

3 

P3 

Trivalent 

Triangular 

^<cP<A 

8 ^ 

P4t 

Tetravalent 

Tetrahedral 

\<cfi<A 

2{A-(f) 

5 

P4s 

Tetravalent 

Square 

l<(fi<A 

(4-d2) 

P6 

Hexavalent 

Cube 

2 < < 4 

(A-d^) 

2 


volume, surface area, and circumference: 

Vd = T2(^ + 4)(2-c?)^ , 

Sd = n{2-d) , 

Crf = iz\/A — cP-. 

B. Overlap Algorithm for Dimpled Spheres 

Standard methods (see, e.g., (29l) for determining the over¬ 
lap between anisotropic bodies rely on convexity. Previous 
theoretical work on concave particles in 3D (TSj [161 EHl ED 
used modifications of an overlap algorithm introduced by He 
and Siders (^ that required the radius of the valence sphere 
to be smaller or equal to the parent sphere and the particles 
were restricted to one valence sphere, {no = 1). For the multi- 
dimpled particles we study here, we use a recently developed 
overlap algorithm (^ . The algorithm (^ is a recursive al¬ 
gorithm based on the Cayley-Menger volume determinant and 
similar dihedral angle determinants. This overlap algorithm 
considers all topological types of overlap among the spheres 
and calculates the mutual union and the mutual intersection 
among them. 

C. Dense Packing 

We vary the square of the distance between the centers of 
the central and the valence spheres in increments of 0.01 
for each of the particles, as shown in Table|T| The domain for d 
varies between particles, so we normalize by linearly mapping 
d^ to the “dimpling amount” /, where / = 0 (/ = 1) when d^ 
is maximum (minimum). 

We construct dense packings of multi-dimpled particles nu¬ 
merically following a standard protocol used in several other 
recent studies (34^33. The densest packing structure is found 
by compressing systems with up to eight particles placed in 
a box with periodic boundary conditions EH EH. Isobaric 
(NPT) simulations, in which the box is allowed to take the 
shape of an arbitrary parellelepiped, are then run with a slowly 
increasing pressure. Each system is compressed 400 times, 
and the densest packing overall is recorded. 

As has been done previously for particles of other shape, 
(351 [39l where tractable, we use numerical putative dens¬ 
est packings as ansatzes for analytical packing constructions, 
which can be used to compute exact dense packing curves. 
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We identify crystal structures by replacing each particle by 
a point at its centroid as in several previous works (TTl |36]- 

Ellllol. 


D. Self Assembly 


We study the self assembly of multi-dimpled particles with 
MC simulations using techniques from several other recent 
studies EElElHia. We use isochoric (NVT) MC simula¬ 
tions with periodic boundary conditions to study systems of 
N particles, with N ranging from 512 to 2064. These system 
sizes are typical of recent work examining entropy driven self 
assembly of hard shapes (e.g. CI][3l|36l|32l) Each system 
is initially equilibrated as a dense, metastable fluid for 10^ 
MC sweeps. Isochoric (NVT) simulations are then performed 
where the packing fraction is slowly increased until it reaches 
a target value. The system is then allowed to equilibrate at 
this state point for ^ 5 x 10^ MC sweeps. Equilibration times 
were chosen based on typical equilibration times that have 
been observed in other studies of hard shapes |[34j [36l, as 
their assembly originates in similar entropic forces. O A 
few experiments were also performed to ^ 10^^ MC sweeps 
to verify that despite their concave shape, no nucleation or 
crystal growth is found to occur beyond typical equilibration 
times for hard convex particles. Densities are varied between 
0.25 < 0 < 0*, where (j)* is the densest known packing frac¬ 
tion for that particle. A full MC sweep for a system of N par¬ 
ticles consists of A -f 1 trial moves including arbitrary transla¬ 
tion, rotation, and box shearing moves. Maximum step sizes 
are periodically adjusted to keep the acceptance probabilities 
at 30%. To check each self-assembled structure for equilibra¬ 
tion, we allow the system to evolve over time scales that are at 
least an order of magnitude larger than the time scale for nu¬ 
cleation and growth, and longer than system auto-correlation 
times. All flndings are verifled by running independent simu¬ 
lations with different initializations and in different box sizes. 
The assemblies reported are those that form at the minimum 
density. 


III. RESULTS AND DISCUSSION 


A. Densest Packings 


Eollowing (3611391 im, we study continuous shape defor¬ 
mations. We compute putative densest packings for each fam¬ 
ily of particles at 100 different dimpling amounts / G [0,1]; 
see Eig. To verify the putative densest packings predicted 
numerically via MC simulations in this section, we perform 
analytic packing calculations for our particles (see appendix 
We And that our numerical calculations match very well 
with the analytic flndings, with an error less than 0.001%. 
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FIG. 2. Numerical calculations of the density (j) of the densest 
obtained packing for dimpled particles as a function of dimpling 
amount /. The curves for different particles are shifted along the 
y-axis for clarity. Packing fraction is maximum at a critical dimpling 
amount / = f^p when the dominating features switch from convex to 
concave regions of the particles. For small /, particles shear from a 
thermodynamic preference for FCC packing in the hard sphere limit 
to a BCT packing with the introduction of dimpling. This is not 
smoothly captured in numerical calculations, but we show that the 
transition is smooth through analytic packing calculations. (Right In¬ 
sets) For P4t and P6 particles, a second transition at larger / is seen 
from /3-tin to FCC and sheared cubic to simple cubic respectively. 
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1. FCC-BCT Transition 


In all cases, for small / ~ 0, particles pack most densely 
into FCC lattices, like hard spheres ( 23 . For slightly larger 
values / < all particle types pack like soft spheres ITt]! 
into a BCT structure. This transition from FCC to BCT pack¬ 
ing structures stems from a new contact between particles 
upon dimpling, dimple-sphere contact, in addition to sphere- 
sphere contact. The dimple-sphere contact reduces the contact 
distance between neighboring particles. This results in a BCT 
structure with a = j3 7 ^ 7 , where a, j3, and 7 are lattice dis¬ 
tances. The transition from FCC to BCT structures occurs 
almost instantaneously with dimpling, and is captured in our 
analytic densest packing calculations. We represent this with 
a single continuous curve in our analytic calculations in Figs. 
Tg[T4 However, in numerical putative densest packing cal¬ 
culations, the particles pack in the FCC structure until a large 
enough dimpling amount due to local frustrations. At these 
values of /, the dimples become increasingly filled by adja¬ 
cent particles with decreasing system volume, allowing denser 
packings without the need for global structural rearrangement. 


2. Critical Dimpling Amount 


We observe a maximum in the densest packing curve of 
each family at a different critical dimpling amount f^p. In 
each family, below the critical dimpling amount / < f^p par¬ 
ticle dimples make contact with the convex part of adjacent 
particles. At / > f^p, neighboring particles interlock within 
each other resulting in a crossover to a structure determined 
by the number and arrangement of dimples on individual par¬ 
ticles. Coincident with the change in the contact between par¬ 
ticles is a crossover in global structure from BCT structures 
to structures that depend on the number and arrangement of 
the concave particle features. Despite an apparent geomet¬ 
ric similarity in the nature of the crossover, i.e. in each case 
the crossover is from BCT to a symmetry dependent struc¬ 
ture, the point at which the crossover occurs, / = f^p, varies 
by approximately a factor of two across particle types. This 
variation suggests that the dimpling amount / is not a well- 
chosen parameter for understanding the origin of the structural 
crossover. Instead, we find that the quantity 


Q 


npCd 

Ind ’ 


( 1 ) 


where npi^ the number of dimples, cj is the dimple circumfer¬ 
ence, and d is the dimple depth, is in the range Cj = 1.23±0.06 
at the crossover across particle types. We note that the depen¬ 
dence of this quantity on as well as cj and J, suggests that 
the crossover is driven by changes in the geometry of the en¬ 
tire particle, and contrasts with what we find for self-assembly 
in section HIT^ below. 

We find that just above f^p, P2 particles pack into a par¬ 
quet structure (BCC without shearing); P3 particles pack into 
an AB-stacked hexagonal packing (hexagonal close packing- 
HCP); P4t particles pack into a j3-tin lattice (sheared dia¬ 
mond); P4s particles pack into a tetragonal lattice with rota¬ 



(a) P2 (b) P3 (c) P4t 



(d) P4s (e) P6 


FIG. 3. Configuration of neighboring particles at the critical dim¬ 
pling amount for densest packing. We observe that the densest pack¬ 
ing configuration is determined by a competition of the circumfer¬ 
ence and depth of the dimples. 


tional parquet symmetry; and P6 particles pack into a simple 
cubic lattice. A further increase in / beyond f^p introduces 
larger voids between particles and a reduction in packing frac¬ 
tion. As / increases beyond f^p, dimple volume continues 
to increase, resulting in a smaller volume of the particle in 
the unit cell structure. This change results in a reduction in 
the packing density of these structures. At sufficiently large 
values of /, P3, P4t, and P6 particles pack less densely than 
spheres. 

For P2, P3, and P6 particles the packing fraction changes 
smoothly as a function of / about its maximum at f^p\ 
whereas for P4t and P4s particles there is a cusp. Both types 
of behavior have been observed previously in dense packings 
of continuously varying shapes ESI EH ED. For P2, P3, and 
P6 particles, the smooth behavior occurs because structures 
on either side of f^p are related by a continuous shear. For 
P4t particles the cusp is the result of a crossover from BCC 
to j3-tin, and for P4s it is the result of a change from BCC to 
simple cubic with a rotational parquet symmetry in the simple 
cubic lattice, which cannot be obtained by a simple shear from 
the BCC structure. 

For P4t particles at high dimpling amounts, / = 0.7374, 
we find another transition from j3-tin to FCC, shown in the 
first right inset of Fig. At these large dimpling amounts, 
the densest packing structure arises from a competition be¬ 
tween the parallel and anti-parallel alignment of the dimples. 
With increased dimple size, neighboring particles find more 
room to rotate while they are interlocked. This results in a 
denser packing when the dimples align parallel to each other 
instead, resulting in a transition from j8-tin, where dimples 
exhibit anti-parallel alignment, to the FCC structure. 

Similar to P4t particles, we find another transition for P6 
particles at / = 0.6510, shown in the second right inset of Fig. 
1^ The particles transition from a sheared cubic arrangement 
to a simple cubic arrangement. This transition is observed 
in hexavalent particles because of their ability to shear along 
the 110 lattice vector direction, while neighboring particles 
remain interlocked in the same configuration. 
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FIG. 4. Packing fraction (j) vs the total dimple circumference-depth 
ratio Q. The packing fraction reaches a critical maximum at = 
1.234±0.060. 

B. Self-Assembly 

We also performed MC simulations at finite pressures to 
study assembly behavior away from the dense packing limit. 


7. Symmetry Independent Behavior near f = 0.- Hard to Soft 
Transition 

In all cases, at / = 0 we observe the assembly of FCC lat¬ 
tices as we would expect for hard spheres. In all cases, as 
dimpling increases from / = 0 we observe a crossover from 
FCC to sheared BCT structures. The onset of the BCT struc¬ 
ture occurs almost instantaneously as a result of continuous 
shearing from FCC into BCC structures. l(42l The FCC-BCT 
transition is followed by a transition from BCT to BCC at dim¬ 
ple volume vj ^ 0.03, which can, again, occur through simple 
shearing of the lattice. ll42l The self-assembly of BCC struc¬ 
tures from nearly spherical hard particles has been observed 
previously (see, e.g., |[II1|371). Here, by means of the radial 
distribution function, g(r), (Fig. shows an example for P4t 
particles) we see that in the BCC phase, neighboring particles 
penetrate into the dimples, but potential of mean force and 
torque (PMFT) calculations (Fig. [7^,c,e and Fig. [^,c) show 
that, overall, the distribution of nearest neighbors remains rel¬ 
atively isotropic despite the presence of dimples. Other work 
(2811431 has shown that soft spheres assemble BCC and other 
non-close packed lattices. In light of this prior work, our ra¬ 
dial distribution function and PMFT results suggest that the 
BCC phases we observe here for nearly spherical particles, 
as well as those observed elsewhere (e.g. imiiTi), arise be- 


Radial Distribution Function 



Radial Distance 


FIG. 5. The radial distribution function of tetrahedrally tetravalent 
particles (P4t) at infinitesimally small dimple volume (above) and 
small dimple volume (below). We report a transition from a face- 
centered cubic to a body-centered cubic crystal structure as volume 
of the dimple increases. 

cause nearly spherical hard particles behave effectively as soft 
spheres, regardless of their symmetry. 


2. Symmetry Dependent Behavior at Large Dimpling 

Above a critical value f>f^ for each family, the phase be¬ 
havior is dependent on the arrangement of dimples, as shown 
in Fig.[^(for futher details, see appendix Figs. p^p^ . P2 par¬ 
ticles continue to self-assemble a BCC lattice, but the lattice 
takes on a parquet form in which particle orientations are or¬ 
dered (see section [niB 3 1 ). For P3 and P4t we do not observe 
spontaneous assembly of an ordered structure from a homoge¬ 
neous fluid at any density on the time scale of our simulations. 
For P3, we find the expected triangular lattice is stable against 
melting to a packing fraction of 0.35 at dimpling amounts 
greater than 0.45. For P4t, if seeded, particles assemble a di¬ 
amond lattice as found previously for tetrahedrally patterned 
enthalpically patchy particles (401 , but unlike truncated tetra- 
hedra (^ or tetrahedrally facetted spheres ITTl . which assem¬ 
ble without seeding. The diamond lattice does not melt down 
to a packing fraction of 0.26 at dimpling amounts greater than 
0.60. P4s particles self-assemble a cubic lattice with rotational 
parquet symmetry. P6 particles self-assemble a simple cubic 
lattice. In all cases we observe that the crossover from BCC 
coincides with an emergent valence, or angular specificity, in 
the effective interactions between neighboring particles, on 
the order of a few k^T (see section [TTlB 3| ), as has been ob¬ 
served for convex particles m 
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FIG. 6. Self-assembled crystals for different multi-dimpled 
spheroidal particles vs volume of individual dimples in particles. 
We find that critical volume of an individual dimple determines the 
spherical (v^ < 0.07) or non-spherical self-assembly of these parti¬ 
cles. We also observe that hard spheres (vj = 0) assemble into the 
FCC lattice. With the introduction of an infinitesimally small dim¬ 
ple, the sheared body-centered tetragonal structure is found. Above 
vj > 0.03, we find that all particles assemble into the BCC lat- 
tice.Error bars are within the size of the gold markers. 


3. Potential of Mean Force and Torque 

The PMFT is a measure of the effective directional entropic 
forces (DBFs) that emerge for a pair of crowded particles a. 
The PMFT (Figs. |7j shows how the local particle environ- 
ment of multi-dimpled lock-and-key colloids differs in differ¬ 
ent bulk structures. We computed the PMFT at the lowest den¬ 
sity where the crystal structure is formed, around a constant 
volume with a single particle placed at the center using meth¬ 
ods introduced in our previous works aim. In the case of P3 
and P4t particles, where no crystal structures are formed, we 
use the lowest stable density of the crystal structure obtained 
from melting studies of the densest packing structure. Figs. [7] 

show isosurfaces of the PMFT in three dimensions plotted 
with VisIT El. In both columns, the PMFT is computed 
on the same constant volume. Below the critical dimpling 
amount /^, we observe a nearly isotropic effective entropic 
potential around the particle. At higher dimpling amounts, 
the dimples behave as an entropic patch ifTTIl that results in 
a preferred directionality between pairs of crowded particles, 
similarly to what was observed previously for facetted spheres 

Emu. 

For large dimples, the PMFT (Figs.|^,d,f and|^,d) show 
that a particle A adjacent to a dimpled particle B will coor¬ 
dinate at a dimple of particle B. This is similar to what has 
been observed previously for facetted spheres and polyhedra 
El El. However, in contrast to what has been observed pre¬ 
viously for other shapes, depending on the structure, particle 
A will either preferentially align its convex region with the 
dimple of particle B, or will align its dimple with the dimple 
of particle B. This is a new feature that arises from the parti¬ 
cle concavity. However, even in the case where both patches 
on neighboring particles align, the behavior differs from what 


has been observed previously. In the diamond crystal struc¬ 
ture obtained from tetrahedrally tetravalent (P4t) particles, the 
dimples align facing each other, similarly to what has been 
observed previously Emu. However, because of the concave 
nature of the entropic patches, there is considerable empty 
volume between the particles, resulting in stability at lower 
packing fraction of 0 = 0.32 than is evident from the structure, 
and what has been observed for convex particles ifTTlITbll . 


4. Crossover in Self-Assembly from Symmetry Independent to 
Symmetry Dependent Behavior 

As in dense packing, with increasing dimple size, in each 
family we observe a change in system behavior. As we noted 
above, in all cases below some /^, we observed the self- 
assembly of BCC phases. Above /^, we observed phase be¬ 
havior that depended on particle symmetry. The fact that in all 
cases we observed a loss of rotational entropy (captured in the 
PMFT results in Figs. 0© driven by changing the particle 
geometry, suggests that there may be a common underlying 
geometric origin of this change in behavior. Depending on 
particle symmetry varied between 0.42 and 0.70 (66.6%), 
which suggests that, as in dense packing, / is not a good pa¬ 
rameter for characterizing the origin of the crossover in bulk 
behavior. Instead, we found that the volume of an individual 
dimple on a multi-dimpled particle showed relatively small 
variation (0.0714 < vj < 0.0783, or Avj/vj « 10%) at the 
crossover in bulk behavior across different particle symme¬ 
tries. 

We found that the change in bulk behavior coincides with 
the binding of particles at dimples, and that this coincides with 
an approximately constant volume of individual dimples. To 
understand why this might be the case, we note that Ref. ||5l 
argued that the phase behavior of several systems of convex 
particles was driven by local packing considerations, which 
in example systems could be reduced to the consideration of 
simple, pairwise packing. Our results suggest the argument 
advanced in Ref. |[5l, i.e. that pairwise packing considerations 
are sufficient to understand the phase behavior of certain con¬ 
vex particles, extends to multi-dimpled sheres. We hypothe¬ 
size why this is the case by considering what happens when a 
particle A binds entropically at a dimple of particle B. When 
particle A binds at the dimple of particle B, particle B loses ro¬ 
tational entropy on the order of ^ In O (where Cl is the number 
of distinct particle orientation^. This entropy loss is com¬ 
pensated by a gain in entropy for the rest of the system. It 
was argued in Ref. O that this should be determined by the 
stress tensor, which we can estimate as typically on the order 
of k^T jvp where Vp is the particle volume. This suggests that 
the additional entropy gained by the system is on the order of 
^BVj/vp, where, as particle pairs bind, we assume that the rest 
of the system gains free volume proportional to the dimple 
volume. There are key differences between systems: the pres¬ 
sure at the onset of ordering, the number of non-degenerate 


^ In technical terms, D is given by the volume of the relevant rotation group. 
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FIG. 7. Potential of mean force and torque (PMFT) calculations fS) show that below the critical dimpling amount (left) an isotropic 
potential exists. At higher dimpling amount (right), an attractive free energy well is presented in the volume of the dimple, giving rise to 
different crystal structures, (a, b) Bivalent particle showing (a) BCC and (b) parquet potentials, (c, d) Trivalent particle showing (c) BCC and 
(d) triangular sheet potentials, (e, f) Tetrahedrally tetravalent particle showing (e) BCC and (f) tetragonal diamond potentials. 
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b Simple Cubic 



FIG. 8. Potential of mean force and torque (PMFT) calculations fS) show that below the critical dimpling amount (left) an isotropic 
potential exists, (a, b) Square Tetravalent particle showing (a) BCC and (b) simple cubic potentials, (c, d) Hexavalent particle showing (c) 
BCC and (d) simple cubic potentials. As the lattice vector of the cubic lattice formed is small in this case, we see second and third neighbor 
peaks in the cubic potential. 


rotations (which is determined by the discrete symmetry of 
the particle), and the orientation of particles coordinated at 
dimples, it would be surprising to find perfect collapse for any 
geometric masure of ordering. Nevertheless, we find the vari¬ 
ation in critical dimple volume is small. 

We expect the result of dimple volume as the controlling 
feature for anisotropic assembly is specific to multi-dimpled 
lock-and-key colloids, and we do not expect this to be a phe¬ 
nomenon that generalizes to all shapes, e.g. to colloidal el¬ 
lipsoids. Indeed prior work on facetted spheres suggested 
there that facet edges were important, but that was based on 
fewer examples than the present work ifTTl . However, what 
we have found is that the suggestion that the phase behavior 


of anisotropic colloids is controlled by pairwise effective in¬ 
teractions, as argued in Ref. O, extends beyond the convex 
particles studied in a, and suggests that the argument ad¬ 
vanced in that work may be universal, even if the quantitative 
details of how it plays out for particles of different types is 
not. 


IV. CONCLUSION 

We found a structural crossover in packing that corresponds 
to a quantity that depends on all features of a particle shape, 
and one for the local ordering in assembly that depends on in- 
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dividual particle features. This finding is interesting in light 
of several studies showing discrepancies between the packing 
and assembly behavior of hard particles (e.g. |[34j[36j[37l[45l). 
Recent work O has argued that the discrepancies between 
packing and assembly arise because assembly is driven by di¬ 
rectional entropic forces that arise from balancing shape en¬ 
tropy, which in turn drives local packing whereas packing at 
infinite pressure is a global phenomenon in which entropy is 
irrelevant. Our finding, that the local ordering at low pres¬ 
sures where we observe assembly (to within 10%) controlled 
the geometry of individual particle features, whereas packing 
behavior is controlled (also to within 10%) by the geometry of 
all particle features, further supports the argument advanced in 

Q. 

Interestingly, multi-dimpled particles self-assemble porous 
crystal structures at packing fractions of 40 — 50% (and un¬ 
der certain conditions can remain at least metastable to much 
lower packing fractions), lower than that observed in convex 
systems (typically > 50% (IT] |36l |37l |45l (461). It has been 
proposed recently that entropy can stabilize open lattices in 
systems with both enthalpic and entropic contributions to the 
free energy 1471 . Here, we found that entropy alone can sta¬ 
bilize porous structures. Intuition might suggest that the ab¬ 
sence of strong enthalpic interactions in the present system 
may also result in fewer kinetic traps, and in several cases we 
found robust bulk 3D assembly. However, the kinetics of lock- 
and-key colloids are subtle: it has been shown that in the pres¬ 
ence of depletants single-dimpled particles first make convex- 
to-convex contact before particles rearrange into convex-to- 
concave contact 1^ . The existence of this complex kinetic 
pathway may explain why we did not observe the spontaneous 
global ordering of P3 or P4t particles here, despite those par¬ 
ticles having DBFs in the range that causes other systems to 
self-assemble. Despite these exceptions, porous structures of 
the type we report here find applications in advanced catalysis 
and medical diagnosis 1^ . The application of recently devel¬ 
oped techniques for optimizing shapes for structures ifT^ to 
multi-dimpled particles might be useful for producing shapes 
that are less likely to get kinetically trapped. 

Finally, we have described a set of systems that exhibit 
complex emergent, entropy-driven phase behavior that can be 
traced to a microscopic quantities that are approximately in¬ 
variant under changes of particle symmetry. Because emer¬ 
gent behaviors are inherently difficult to ascribe to micro¬ 
scopic system details, we believe this example in experimen¬ 
tally realizable systems of colloids OUEll will provide a use¬ 
ful experimental setting for manipulating emergent behavior. 
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Appendix A: Geometric Characteristics of Multi-Dimpled 
Particles 

We calculate geometric characteristics of the particles at 
various dimpling amounts, as shown in Figure 

1. Dimple Volume 

The volume of a single dimple, vj, is the union of the vol¬ 
ume of the two intersecting spherical caps that form the dim¬ 
ple. The volume of spherical cap of radius R and height h is 
given by 

V{R,h) = Inh^i^R-h). 

For a single dimple formed at the intersection of central and 
valence spheres of radii R+,R_, its volume is given by 

~ (^+ + ~ X 

(J^ + 2JR+ + 2dR- + 6R+R_ - ?>r\ - 3rI) (A1) 
. For R+ =R- = 1, 

Vd=^(d + 4)(2-df. (A2) 

where d is the distance between the centers of the central 
sphere and the valence sphere. The volume of the remaining 
particle Vsp is given by 

Vsp = vps - {no X Vd). (A3) 

where vps is the volume of the original central sphere, vj is the 
volume of a single dimple and is the number of dimples in 
the particle. This gives 

v^ = ^nRl-nD-^^{d + A){2-df. (A4) 

We normalize the total dimple volume for a particle by 


^sp 

For normalized single dimple volume, we set /t/) = 1. 

2. Dimple Surface Area 

Similar to the volume of the dimples, we calculate the sur¬ 
face area of a single dimple and normalize the total surface 
area of the dimples by the surface area of the particle. The 
surface area of a single dimple is calculated as the surface 
area of a spherical cap formed between the central and valence 
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(a) Circumference (b) Surface Area (c) Volume 


Appendix B: Critical Dimpling Amount 


TABLE II. Crossover Values for Putative Densest Packings and Ther¬ 
modynamic Assemblies 


FIG. 9. The different shape features of multi-dimpled spheres con¬ 
sidered in this study. 


spheres. For spheres of radii the surface area is given 

by 

Sd='^{R+-R-+d){R++R--d). (A6) 

For=/?_ = !, 

Sd = 7t{2-d). (A7) 

The surface area of the particle, Ssp, is given by 

Ssp = sps - {no xsd). (A8) 

where sps is the surface area of the original central sphere. 
This gives 

s^p = 4nR\ -hdX 7r{2 - d ). (A9) 

The self-assembled non-spherical crystal structure is rather 
dependent on the number and arrangement of dimples in the 
particle. We normalize the total dimple surface area Sd for a 
particle as 




np xsd 

^sp 


(AlO) 


3. Dimple Circumference 


Similarly, the circumference of the dimples cj is calculated 
from the radius of the circle of intersection between the central 
and valence spheres. For spheres of radii /?+,/?_, the circum¬ 
ference is given by 


Cd — 2 \/ (~^+ + (^+ “ X 

2{R+ dR+-d) (/?+ dR-dd) (All) 

. For/?+ =/?_ = !, 

Cd = n\/‘\-d^. (A12) 


This circumference is normalized by the distance between the 
central and valence spheres. To reduce the numerical ratio, 
we use a weight 27r and call this normalization factor Csp. This 
factor is given by 


Csp = 2nd. (A 13) 

The total circumference of the dimples Q is thus normal¬ 
ized by the circumference of a circle describing the distance 
between the parent and valence sphere in the particle. 


Q 


np X Cd 
Csp 


(AI4) 


Particle Q 


P2 

0.73 

1.21 0.42 0.0742 

P3 

0.53 

1.22 0.48 0.0714 

P4t 

0.35 

1.20 0.53 0.0782 

P4s 

0.66 

1.32 0.73 0.0769 

P6 

0.30 

1.26 0.70 0.0783 


The critical dimpling amount for the putative densest pack¬ 
ing and self-assembly of each type of multi-dimpled 
particle is the crossover point in the corresponding particle 
behavior. At these crossover values, different geometric char¬ 
acteristics of each particle type have been computed as shown 
in Table HD 

To understand the packing of multi-dimpled spheroidal par¬ 
ticles, we determine the primary feature in the shape of these 
particles that affects their behavior. We investigate the value 
of compared to the total volume Vj, surface area Sd and 
circumference Q of all dimples on a single particle as d is 
varied. We normalize each of these three parameters. The 
volume and surface area are normalized by the remaining re¬ 
spective quantity in a single particle, while the circumference 
is normalized by the circumference of a circle with radius 
equal to the distance between the central and valence spheres, 
a linear function of the depth of the dimple. We compute the 
total normalized dimple volume, surface area and circumfer¬ 
ence to depth ratio in each family of particle at We find 
that the total circumference to depth ratio Q of the dimples 
is 1.23±0.06 at f^p across all five particle families; see Fig. 
1^ In contrast, no dimple volume or surface area showed no 
similar collapse. 

We also observed that in self-assembly, all particle types 
undergo a structural transition from BCC in which parti¬ 
cles have nearly isotropic pairwise interactions, into a phases 
where particle symmetry is important. To quantify when the 
concave features are sufficient to induce a change in the pre¬ 
ferred structure, we compare particle shape characteristics at 
the critical dimpling amount for self assembly /^. We find 
that, across particle types, vj shows little variation at /^, 
which suggests that the volume of individual dimples controls 
assembly behavior, independent of particle type. 


TABLE III. Geometric and Numerical Calculations of Particle Inter¬ 
locking 


Particle Description f^p fj^p d^ 

P2 Bivalent 0.73 0.86 2.0 

P3 Trivalent 0.53 0.70 2.14 
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Appendix C: Analytic Calculations of Packings 


To verify the putative densest packing predicted computa¬ 
tionally via MC simulations for the particles studied, analytic 
calculations of the packings of these particles were performed. 
In Figures [^-[^ we plot analytic calculations of the pack¬ 
ing curve and numerically predicted packings for one and two 
particles of all particle types. We show that our numerical cal¬ 
culations match the analytic findings to within 0 . 001 % when 
/ > 0.10. Additionally, we perform numerical and analytic 
packing calculations for single dimpled particles (see Figure 
[Tg. We find that single dimpled particles do not follow the 
pattern shown by other particles in this study due to the lack 
of discrete symmetry in their shape. 


We further numerically calculate maximum packing den¬ 
sities of all particle types with up to eight particles. How¬ 
ever, due to the complexity of the analytic packing problem 
for number of particles ^ > 3, we do analytic calculations only 
for (some) one- and two-particle packings {n= 1 , 2 ) for dif¬ 
ferent particle types. If solvable, these analytic curves and 
intersection equations are shown in Figures [^-[T^ For these 
results it is useful to work in terms of a linearly shifted param¬ 
eter F = 0.5J, where d is the distance parameter, rather than 
the dimpling parameter / we use in the main text. We denote 
the volume of a particle by U, the volume of a unit cell of the 
lattice by V, the packing fraction by (p, and different packing 
regions by T. We denote the vectors describing a unit cell of 
the packing as a, j3, and 7 , and the volume is given by the 
standard formula V = a • (j3 x 7 ). Images are colored so that 
the color on the left hand side of the image matches the color 
used for the particle in the rest of the text, and then proceeds 
through the visual spectrum over the range of geometrically 
allowed particle parameters starting toward the red end of the 
spectrum for single particle packings, and the blue end of the 
spectrum for two particle packings. 


The relationship between the thickness parameter F and 
dimpling parameter / for each particle type is included in the 
bottom left side of the graph. Analytic calculations of the di¬ 
mensions of the packing box are shown in the graph. The an¬ 
alytic packing density curve is drawn as a line in these figures 
against the thickness parameter F, in regions where it can be 
calculated. In order to compare with putative density calcula¬ 
tions, these lines are superimposed on markers denoting nu¬ 
merical calculations. Various regions in these analytic pack¬ 
ing curves are explicitly marked through vertical lines at dif¬ 
ferent values of F, the thickness parameter. In regions where 
the analytic curves can be calculated, we find that they match 
our numerical calculations, with an error less than 0 . 001 % 
when / > 0.10. We include the single dimpled particle PI, 
since the same lattice is shared between particles PI and 
P2. Similarly, the same lattice is shared among particles 
P2, P4s and P 6 . 



FIG. 10. Analytic and numerical packing curves showing a perfect fit 
for bivalent particle (P2) for one and two particles in a box. Analytic 
packing curve for one particle is shown in red, two particles in black. 
The dark markers denote numerical densest packing calculations for 
two particles in a box, and light for one particle. Blue dashed lines 
show the different packing regions found in numerical calculations. 


P3 - Analytic and Numerical Packing Curves 



Dimpling Amount / 


FIG. 11. Analytic and numerical packing curves showing a perfect 
fit for trivalent particle (P3) for one and two particles in a box. The 
curves overlap as one and two particles pack in the same manner. 
Green dashed lines show the different packing regions found in nu¬ 
merical calculations. 
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P4t - Analytic and Numerical Packing Curves 

1 . 0 ^. 



Dimpling Amount / 


FIG. 12. Analytic and numerical packing curves showing a perfect fit 
in a single region between black dashed lines for tetrahedral tetrava- 
lent particle (P4t) of one and two particles in a box. Analytic packing 
equations could not be calculated at all dimpling amounts by analytic 
means. Red dashed lines show the different packing regions found in 
numerical calculations, in addition to those found in analytic calcu¬ 
lations. 


P4s - Analytic and Numerical Packing Curves 



Dimpling Amount / 


FIG. 13. Analytic and numerical packing curves of a square tetrava- 
lent particle (P4s) showing a perfect fit in the entire domain for one 
particle in a box and between black dashed lines in two regions for 
two particles in a box. Analytic packing curves could not be calcu¬ 
lated at all dimpling amounts for two particles by analytic means. 
Magenta dashed lines show the different packing regions found in 
numerical calculations, in addition to those found in analytic calcu¬ 
lations. 


Appendix D: Self-Assembly 
1. Transition Density 


The particles studied in this work assemble at approxi¬ 
mately the same conditions. For the sake of completeness, 
we show the transition density phase diagrams for thermo¬ 


dynamic self-assembly of each particle in Figures 18-22 


In these plots, we also observe the reported behavior of a 
crossover in the assembly of particles beyond a constant criti¬ 
cal volume of a single dimple in the top axis. 

These crystal structures are formed by entropic interactions 
between the convex and concave regions of two adjacent par¬ 
ticles. We find that the emergent valence is dependent on 
the geometry of the dimples (see section [TTl B 3| for details). 
The number of directions for the alignment of the particles 
depends on the crystal lattice, for example, two in the par¬ 
quet structure. It should also be noted that the particles have 
rotational symmetries, which align with the symmetry of the 
crystals and thus do not break the symmetry of the crystals. In 
the case of the tetrahedrally tetravalent particle, particles have 
two orientations in the diamond crystal structure formed. 


The self-assembly of multi-dimpled spheres into crystal 
structures occurs at packing fractions of 40 — 50%, much 
lower than that observed in entropically driven assemblies of 
convex polyhedra, typically > 50%. These assemblies are 
thus candidate structures for entropically assembled porous 


structures, and can find applications such as advanced cataly¬ 
sis. By designing porous structures using entropy alone, their 
assembly pathways have lesser kinetic traps than their en- 
thalpic counterparts, giving rise to structures that can be easily 
assembled and reconfigured. 


2. Dimpling Amount and Volume of a Dimple 

As shown in Figures the phase diagrams do not 

collapse at a constant dimpling amount /. However the vol¬ 
ume of a single dimple is constant at the crossover across all 
phase diagrams. Here we include calculations that show the 
relationship between volume of a single dimple and the dim¬ 
pling amount, in order to show that these two quantities are 
not linear functions of each other. 

For bivalent particle, P2, with two dimples, we have: 


/= 


4-1 


(Dl) 


d^=4-3f. 


where, / is a linear mapping of cfi, the distance between the 
valence and central spheres. 
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P6 - Analytic and Numerical Packing Curves 

1.0 


0.9 


.2 0.8 

I 

pL, 

I 0.7 

c?J 


0.6 


0.5 

0.0 0.2 0.4 0.6 0.8 1.0 

Dimpling Amount / 

FIG. 14. Analytic and numerical packing curves of a hexavalent par¬ 
ticle (P6) showing a perfect fit between red dashed lines in two re¬ 
gions for one particle in a box and between black dashed lines in 
one region for two particles in a box. Analytic packing curves could 
not be calculated at all dimpling amounts for two particles by ana¬ 
lytic means. Orange dashed lines show the different packing regions 
found in numerical calculations, in addition to those found in analytic 
calculations. 



Next, we have: 


{d + A){2-dY 

%-{d + A){2-dY ■ 


(D2) 


where, Vj is the normalized single dimple volume, equivalent 
to the ratio of the volume of a single dimple to the volume of 
the remaining particle. 


In a similar fashion, we understand that the relationship be¬ 
tween the volume of a single dimple and the dimpling amount 
is different for each particle type and is a function of the dim¬ 
ple geometry in a particle. 
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FIG. 15. Analytic packing density ^ < 0 < 1 and packing curves as a function of thickness parameter 0 < F < 1 (right) for single & double 
lattice packings of particles PI & P2. 
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FIG. 16. Analytic packing density ^ < 0 < 1 and packing curves as a function of thickness parameter 0 < F < 1 {right) for single & double 
lattice packings of particles P3 & P4t. 
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FIG. 17. Analytic packing density ^ < (j) < I and packing curves as a function of thickness parameter 0 < F < 1 (right) for single & double 
lattice packings of particles P4s & P6. 
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FIG. 18. Minimum packing fraction at which bivalent particles (P2) 
assemble crystal structures. The structures assembled vary with dim¬ 
pling amount /, from FCC for hard sphere to BCT, BCC and a par¬ 
quet structure. Blue dashed lines show the different structural regions 
found in our simulations. 


FIG. 20. Minimum packing fraction at which tetrahedral tetravalent 
particles (P4t) assemble crystal structures. The structures assembled 
vary with dimpling amount /, from FCC for hard sphere to BCT, 
BCC and diamond. Red dashed lines show the different structural 
regions found in our simulations. 
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FIG. 19. Minimum packing fraction at which trivalent particles (P3) 
assemble crystal structures. The structures assembled vary with dim¬ 
pling amount /, from FCC for hard sphere to BCT, BCC and then 
form glass. Green dashed lines show the different structural regions 
found in our simulations. 


FIG. 21. Minimum packing fraction at which square tetravalent par¬ 
ticles (P4s) assemble crystal structures. The structures assembled 
vary with dimpling amount /, from FCC for hard sphere to BCT, 
BCC and simple cubic. Purple dashed lines show the different struc¬ 
tural regions found in our simulations. 
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FIG. 22. Minimum packing fraction at which hexavalent particles 
(P6) assemble crystal structures. The structures assembled vary with 
dimpling amount /, from FCC for hard sphere to BCT, BCC and sim¬ 
ple cubic. Yellow dashed lines show the different structural regions 
found in our simulations. 
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